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ABSTRACT 
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• Context. A number of binary systems present evidence of enhanced activity around periastron passage, suggesting a connection between tidal 
I interactions and these periastron effects. 

O Aims. The aim of this investigation is to study the time-dependent response of a star's surface as it is perturbed by a binary companion. Here 

we focus on the tidal shear energy dissipation. 
^ ^ Methods. We derive a mathematical expression for computing the rate of dissipation, E, of the kinetic energy by the viscous flows that are 
I I . driven by tidal interactions on the surface layer of a binary star. The method is tested by comparing the results from a grid of model calculations 

' with the analytical predictions of Hut (1981) and the synchronization timescales of Zahn (1977, 2008). 
T-H Results. Our results for the dependence of the average (over orbital cycle) energy dissipation, Eay,,, on orbital separation are consistent 
^ with those of Hut (1981) for model binaries with an orbital separation at periastron rp^JRi > 8, where Ri is the stellar radius. The 
model also reproduces the predicted pseudo-synchronization angular velocity for moderate eccentricities (e <0.3). In addition, for circular 
orbits our approach yields the same scaling of synchronization timescales with orbital separation as given by Zahn (1977, 2008) for 
^0 convective envelopes. The computations give the distribution of E over the stellar surface, and show that it is generally concentrated at the 
equatorial latitude, with maxima generally located around four clearly defined longitudes, corresponding to the fastest azimuthal velocity 
perturbations. Maximum amplitudes occur around periastron passage or slightly thereafter for supersynchronously rotating stars. In very 
eccentric binaries, the distribution of E over the surface changes significantly as a function of orbital phase, with small spatial structures 
' appearing after periastron. An exploratory calculation for a highly eccentric binary system with parameters similar to those of 6 Sco 
" _ I . (e=0.94, ^=3944. 7 d) indicates that changes by ~5 orders of magnitude over the 82 days before periastron, suggesting that the sudden and 
^ ' large amplitude variations in surface properties around periastron may, indeed, contribute toward the activity observed around this orbital phase. 
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1. Introduction 

A number of binary systems present evidence of enhanced ac- 
tivity around periastron passage. Among these, 77 Car and the 
Wolf-Rayet systems WR 140 and WR 125 are the most ex- 
treme and best documented examples of periodic brightening 
at X-ray, visual and IR wavebands associated with periastron 
passages. Recently, van Genderen & Sterken (2007) suggested 
that these periastron events may have the same physical cause 
as the milder "periastron effects" exhibited by many renowned 
eccentric binaries in which small enhancement (Am,, ~0.01- 



0.03"""*) in the visual brightness of the system around pe- 
riastron passage are observed. They suggest that the funda- 
mental cause of the effects may reside in the enhanced tidal 
force that is present during periastron passage. In addition to 
brightness enhancements, binary interactions are frequently in- 
voked to explain certain mass-ejection phenomena. For ex- 
ample, Koenigsberger, Moreno & Cervantes (2002) raised the 
question of whether tidal forces play a role in changing the 
wind structure in the massive Small Magellanic Cloud binary 
system HD 5980, and suggested a possible link between tidal 
forces and the instability producing the 1994 eiTiptive event. 
Bonacic Marinovic et al. (2008) proposed a tidally-enhanced 
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model of mass-loss from AGBs. Observational data supporting 
the idea of mass-loss events associated with periastron passage 
exist for the highly eccentric (e ~0.94) binary 6 Sco (Bedding 
1993; Miroshnichenko et al. 2001; Tango et al. 2009). And 
Millour et al. (2009) suggested that the dusty circumbinary 
environment in HD 87643 might be associated with repeated 
close encounters in the long period and possible highly eccen- 
tric binary. In this paper we explore a mechanism that has the 
potential of providing a theoretical framework for analyzing 
these phenomena. 

Tidal interactions are ubiquitous in binary systems of all 
types. They lead to the deformation of the stellar surface. When 
the system is in equilibriurrQ, the deformation remains con- 
stant. However, if the stellar rotation is not synchronized with 
the orbital motion, the shape of the star is time-dependent. In 
eccentric binaries, the orbital separation changes as the stars 
move from periastron to apastron, thus leading to a time- 
dependent gravitational force. In addition, because the orbital 
angular velocity, Q, is a function of the orbital separation, the 
departure from synchronicity between the stellar rotation an- 
gular velocity, w, and Q is orbital-phase dependent. Hence, ec- 
centric binaries are never in synchronous rotation. The asyn- 
chronicity leads to the appearance of potentially large horizon- 
tal motions on the stellar surface, which are referred to as tidal 
flows. 

The energy that is dissipated because of shearing motions 
is ultimately converted into heat that is deposited in the stellar 
layers. Whether this process leads to detectable observational 
effects is one of the questions that has driven our investiga- 
tion of the detailed tidal interaction effects. In the first paper 
of this series, a very simple model was presented for comput- 
ing the stellar surface motions in a component of a binary sys- 
tem (Moreno & Koenigsberger 1999, hereafter Paper I). In that 
model we computed only the motion of surface elements on 
the equatorial plane of a star. We followed a Lagrangian ap- 
proach in a quasi-hydrodynamic scheme, solving the equations 
of motion of small surface elements, as they respond to gravi- 
tational, gas pressure, viscous, centrifugal, and Coriolis forces. 
The model was applied to the l Ori binary system (Paper I), 
the e Per system (Moreno, Koenigsberger & Toledano 2005, 
hereafter Paper II), and the optical counterpart of the X-ray bi- 
nary 2S01 14+650 (Koenigsberger et al. 2006). In Paper II we 
computed the absorption line profiles obtained with the model, 
taking an ad-hoc extension to stellar polar angles of the motion 
computed on the stellar equator. This extension was also used 
in the analysis of energy dissipation rates and synchronization 
time scales in binary systems with circular orbits (Toledano et 
al. 2007). Here we improve the model by computing explic- 
itly the motion of elements along different parallels that cover 
the stellar surface. This new scheme in the computation of en- 
ergy dissipation was applied in the LBV/WR system HD 5980 
(Koenigsberger & Moreno 2008), and in the study of tidal flows 
in a Virginis (Harrington et al. 2009). 



' A binary star is defined to be in equilibrium (Hut, 1980) when it is 
in a circular orbit (e=0), its rate of equatorial rotation equals the rate 
of orbital motion (ai, ,,, = Q.) and the axes of stellar and orbital angular 
velocity are both perpendicular to the orbital plane. 



This paper is organized as follows: In Section 2 we gener- 
alize the mathematical expression for the calculation of shear 
energy dissipation rates, E, that was presented in Toledano et 
al. (2007); in Sections 3 and 4 we compare the predictions of 
our model for E with the analytical results of Hut (1981) and 
the synchronization timescales with the analytical expression 
of Zahn (1977, 2008), respectively; in Section 5 we explore 
the orbital-phase dependent effects and their possible relation 
to the observational periastron effects; in Section 6 we present 
a discussion, and Section 7 lists the conclusions. 

2. The extended model 

As in the first version of the model, we assume that the main 
stellar body, below the thin surface layer, behaves as a rigid 
body. In the quasi-hydrodynamic method used in Paper I, for 
every surface element a detailed analysis of the positions of 
neighboring elements is needed to assign a local dimension. If 
an initial well ordered arrangement of the elements begins to 
mix freely in the azimuthal and polar directions, the assign- 
ment of a physical dimension to an element and its correspond- 
ing interactions with neighboring elements is not treatable in 
our current scheme. Thus, we make the simplifying assumption 
that the surface elements move only in the radial and azimuthal 
directions; that is, an element always stays in its initial paral- 
lel, and there is no meridional motionjl The energy dissipation 
computation performed in Section |23] requires the motions of 
the surface elements. In this section we give the procedure to 
compute this motion in either of the two stars in the binary sys- 
tem, say star 1, with mass m\ and initial radius R\. Star 2 has 
a mass m2 and instantaneous position r2\ with respect to the 
center of mi, and an orbital angular velocity Q. The motions 
of the surface elements in m\ are computed in a primed non- 
inertial reference frame with its origin at the center of nii, and 
Cartesian axes rotating with the orbital angular velocity £1; the 
x' axis always points to mi- 

In Paper II the total acceleration a' of a surface element, 
measured in the non-inertial frame, is shown to be 



Gniir' 



r - r2i 



\r' 



-Q X (Q X r') - 2il X v' 



\r2\ 



dt 



X r . 



(1) 



with a* the acceleration of a surface element produced by gas 
pressure and viscous forces exerted by the stellar material sur- 
rounding the element, and r', v' the position and velocity of the 
element in the non-inertial frame. A nearly spherical shape is 
assumed for nii throughout its motion. 

The computation starts at time f = with no relative motion 
between the surface elements, so that at this initial time the 
acceleration a*, with value has only the contribution of gas 
pressure. To compute this acceleration we consider a second. 



^ Eqs. 40b and 40c given in Scharlemann (1981) indicate that for a 
colatitude angle 9 >70°, the maximum perturbation in the azimuthal 
direction is more than 15 times greater than the maximum perturbation 
in the polar direction. 
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double-primed, non-inertial reference frame, whose origin is 
also at the center of mi, and has axes rotating with the assumed 
constant angular velocity of the inner rigid region of star 1 . 
This angular velocity is written as - /3oilo, that is, equal to 
a certain fraction ySo of the initial orbital angular velocityjj The 
condition for initial equilibrium Vq = 0, a'^ = on the surface 
of mi is (with r" = r') 



this azimuthal direction, and is twice the distance between 
the center of mass of the element and the boundary of the inner 
stellar region. At f = all the elements in a given parallel have 
the same lengths. 

The acceleration has contributions from gas pressure 
and viscous shear. Below we describe the components of these 
contributions in spherical coordinates. 



' Gm\r' 



„/|3 



+ Gm2 



r - 1-21 



\r'-r2x\ 



k2il 



f=0 



(2) 



The initial position of a surface element is obtained with 
an initial spherical shape of star 1, and its initial velocity is Vg 
= 06o - 1)^0 X r;. 

The contribution of the m2- term in Eq. (|2]i makes a*,, in a 
given parallel have a non-constant azimuthal component. Thus 
the assumed initial spherical shape would not be the proper 
shape consistent with this azimuthal behavior. In our compu- 
tations, at f = we impose only the radial equilibrium result- 
ing from Eq. (|2]i, and allow unbalanced forces in the azimuthal 
direction. We found that with appropriate parameters (e.g. vis- 
cosity) this non-equilibrium initial condition will be followed 
by a transient phase and a steady state (dependent on the orbital 
phase) at later times. 

At f > the acceleration produced by gas pressure is con- 
veniently written in terms of the radial component of a*o. Thus 
the radial component in Eq. (|2]i, which is used in the following 
section, is 



+ Gm2 



+ Gm2 



r'o -r2io sine' cos 



('"o^ + '"210 ~ 2r;)r2i() sin 6' cos i^^,) 



3/2 



sin ff cos (p'q 



'210 



-jSoQor^jSin^e', 



(3) 



with 0' the polar angle of the given parallel and i^,', the initial 
azimuthal angle of the element. 

2.1. The acceleration a-t. 

An element on the stellar surface of mi has lateral surfaces fac- 
ing in the three directions of spherical coordinates, r', (p', 9' in 
the primed non-inertial frame. The x'- axis points to m2 and (p' 
= on the positive side of this axis; the polar axis is z', which 
is also the stellar rotation axis. We will denote with / the paral- 
lel's number and with j the number of the element. The lengths 
of an element in the three directions are /,' , L' , and L' , with 
initial values , I0' , h' ■ The length Ig' is assumed constant 

i jo '^i jo ijo i j 

in time for all elements in a given parallel. In our scheme, for a 
given element at any time is the mean of its azimuthal dis- 
tances to the centers of mass of the two adjacent elements in 



The expression y8,,f,- = a>Q/Q.p„ = 0.02 P ■ 



^ui ^'■per - "-"^ i ^ convenient 

form to compute the synchronicity parameter at periastron. P is in 
days, v,„, in km s"', and Ri in Rq. 



2.1 .1 . Gas pressure 

The gas pressure inside a surface element is pij, and we as- 
sume a polytropic state equation pij - Pijoipijl PijoY , with 
y = 1 H- l/«; n is the polytropic index and p the mass den- 
sity. The pressure on the element exerted by the neighboring 
stellar inner region is taken as /?,„, = Pij/q, < q < 1. We 
consider in particular the half value q = 0.5. Thus the initial 
radial equilibrium is pijok', l^p'. Iq - mij{ai,or')ij, with niij the 
constant mass of the element. At f > the radial pressure ac- 
celeration is (ai,r'\)ij - Pijle'. lip'. l qfHij. This relation combined 
with the polytropic equation, the initial radial equiUbrium, and 
with the constant value of /«' , gives 



(«*;•' 1 )ii 



lr'.\ 



(O-koi')i 



(4) 



The azimuthal gas pressure acceleration on a given ele- 
ment, (fl*;£,'i), /, is computed with the difference of the gas pres- 
sure forces exerted by the two adjacent elements in this az- 
imuthal direction. Then (a*^'i)/j = (Pij-i - Pi,j+\)lr' k'./mij. 
With the polytropic form for each pressure, the corresponding 
initial radial equilibrium conditions on the adjacent elements, 
and the same initial mass density for each surface element (thus 
the same mass of elements in the given parallel, because their 
initial lenghts are the same), this expression is 



'^'ijo 

X {(flw)/,j-i [/r;, ^ -(fl*or')u+i [^;,^,^^;,^,] ^j- 

(5) 



2.1 .2. Viscous force 

The shear part of the kinematic stress tensor is given by (e.g., 
Symon 1971) 



V'v' + (VV)"''--l(V'-v') 



(6) 



where 1 is the unit matrix, trp indicates the transposed, and 77 is 
the coefficient of dynamic viscosity, related with the coefficient 
of kinematic viscosity, v, by 77 = vp, with p the mass density. 
The total force on a given surface element from shear stresses is 
F;, = - ^P,, ■ dS, with the integration over its surface. We pic- 
ture each element with lateral faces in the three directions of 
spherical coordinates, and corresponding lateral areas Ar', A^p^, 
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Agi . Because we do not compute the motion in the polar direc- 
tion, with w' - v^,/r' sin 6', the relevant radial and azimuthal 
components of the viscous force are 

,dcj' 



+ Ao,A 



'^[3'' dr'[r' j 3 dtp' 
1 dV, 



77 1 r — — sm & H — 

' \ dr' r' sm 6' dip' 



1^ 
r' dO' 



,dcj' . , 1 dV^, 

n I r — — sm 6 H — 

'1 dr' r'smff dip' 



+A^'A 



2 , d h'A Adu' 



''I 3'' dr'\r' 3 dip' 



(7) 



In the first term of Eq. (fTOl l we approximate ido)' jdr') ^ 
(u)' - <jL>'^)llr', with w'^ - ySoQo - the angular velocity of 
the inner stellar region as measured in the primed non-inertial 
frame. In the second term we make a similar approximation as 
that used to obtain Eq. ( fTTT i and ignore the factor 4/3. Thus, the 
corresponding accelerations are 



V 1 



V , 



v' + v' - 2v' 



(13) 



(14) 



and finally, the last term in Eq. (fTOl i gives the azimuthal accel- 
eration 



+ Ag,A 



dcj' 



n— — sin0' 



(8) 



VI ■ ■ 

(a*^'4),7 ^ [i^'ilx,] - ^'i) sin - (w)^. - oj'U,) sin Q\_^\ , (15) 



with A meaning the difiference of values of the function inside a 
square parenthesis computed at the corresponding lateral faces 
in the direction given by the lateral area factor. 

In our macroscopic element description, for the evaluation 
of the several terms in this surface integration we consider as 
important the gradients across the lateral faces of an element 
and ignore local gradients along these faces. Thus, given that 
the density drops to zero at the outer face and the radial velocity 
is also zero at the inner rigid region, we approximate the radial 
and azimuthal components of the viscous force as 



■ A^,A 



5v' 



r' sin 6' dtp' 



+ AffA 



r' d9' 



(9) 



with 0)'* , ., o)'* , . the mean angular velocities of the elements 
in the polar direction adjacent to the element ij, and 6'^_^ ■, 9'..^-^ 
the polar angles of the boundaries of the adjacent parallels. 



2.2. Equations of motion 

The acceleration a* of a surface element produced by gas pres- 
sure and viscous forces has a radial component a*,' given by 
the sum of Eqs. (|4|, (fTTI) . and (fT2l) . and the corresponding az- 
imuthal component a*^- is the sum of Eqs. (|5]l, ( fT3l l. ( fT4l l. and 
(fTsT l. In the primed non-inertial reference frame, the pair of ra- 
dial and azimuthal equations of motion of a surface element 
are 



(F,)^, ^ -A 

+AffA 



r\r — — sm Q 



dr' 

doj' . , 

ri smfc* 

' d0' 



4 



doj' 



(10) 



with the first term in Eq. (fTOl l computed at the boundary with 
the inner rigid region. 

The terms in Eq. (|9]l give second and third contributions 
to the radial acceleration, besides that in Eq. Q. Taking the 
coefficient of kinematic viscosity, v, to be the same in all the 
elements, the acceleration corresponding to the first term is ap- 
proximated as in Eq. (1 1) of Paper II 

(a*.2)7-f K,,+<,^_,-2vg, (11) 

f'.i 

likewise, the second term in Eq. (|9]l gives the acceleration 



(12) 



with v'* , v'* the mean radial velocities of the elements in 

' i-ij r i^ij 

the polar direction adjacent to the element ij. 



r' - r2i sin 6' cos(^' 



sin 6' cos <f' 



(r'^ + r^j - 2r'r2\ sin 6' cos ip'^ 



3/2 



+ (fl + ip'fr' sin' ff. 



(16) 



1 



r' sin 6' 



I a*„' - Gmi 



(r'^ + j-ji - 2r'r2i sin9' cosip'^ ' ''21 



simp 



-Q {n + ip')r, 

y/ 



(17) 



All pairs of equations in all surface elements are solved simul- 
taneously, along with the orbital motion of 1112 around mi. The 
values of rai, SI, Cl are obtained from this orbital motion. We 
used a seventh-order Runge-Kutta algorithm (Fehlberg 1968) 
to solve the equations. With some forty parallels covering the 
stellar surface, more than 10'* elements are employed. 
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2.3. Method for E calculation 

The rate of energy dissipation per unit volume can be expressed 
in terms of the tensor in Eq. (|6), and is given by the matrix 
product (e.g.,McQuarrie (1976)) 



-P, : V'v'. 



which can be written as 
Ey - -2/; 



i(V'.v'f -(VV),:(VV), 



(18) 



(19) 



where (V'v')s is the symmetric tensor (V'v')s - jV^' v'+ 
(Vv')"^]. 

The computations in our model show that the gradient in 
azimuthal velocity dominates that in radial velocity. Thus, Eq. 
ST9[ can be approximated as 




-100 100 

Azimuth angle (p'C) 



1 



3 \r' sin 6' d(p' 



dr' 



r' 



dv' 
f_ 

dO' 



tanfl' 



(20) 



Now with oj' - y'^ijr' sin0', it follows that 

v;, 

-TT- - — = r'— sine', 
or' r' or' 



— i- = /— smfl', 
dip' dip' 



, d(x>' . , 
r — — sm 6 , 
36' 



86' tan 6' 
and Eq. ( l20l l reduces to 



\Aldu' 



dr' 



du)' 



sin^e' 



(21) 
(22) 

(23) 
(24) 



For an accretion disk {ff = n/2, and the primed reference frame 
is taken as inertial in this case) with ai' a function only of the 
radial distance r', and v^, small compared with v^,, Eq. ( l24l i 
gives Ev - w'^i^f (Lynden-Bell & Pringle 1974); this result 
was used in equation (14) of Toledano et al. (2007), multiplied 
by the volume of an element. 

The j- surface element on the parallel with polar angle ff. 
has a volume AV// = r'^Ar' Ato'. A6''sin0' - L' U h' , with 

J IJ IJ ' IJ I I if 'Tij ij 

L' = r' . sin 9'Aip' . and la> - r' AO'.. The rate of energy dissipa- 

Y,j !J I I IJ ,j 'J ' ± 

tion in the element is Ejj - Evij^Vij, with Eyij the evaluation 
of Eq. (l24l i at the position of the element. 

As already done to approximate Eq. (fTOt . we take 



a>.. - 



I,-'.. 



{dan _ HKi+i-Ki-i) 



(25) 



(26) 



Fig. 1. Rate of energy dissipation Ejj as a function of ip' in a 
few selected parallels, for the n=1.5 polytropic case hR\ jRi - 
0.03, V = 0.003 Rl day-\p = 8.66x10"^ gr cm-\ orbital period 
f = 10 days (a - 40.62 Rq), and y6o = 1.2. The top curve shows 
results for the stellar equator, and following curves from top to 
bottom correspond to parallels with 0' 81°, 73°, 64°, 55°, 46°, 
38°, 29°, and 20°. Calculations with /3o = 2.0 yield identical 
curves, but the Eij values are higher for by ~1.5 magnitudes. 



Iff Ir'. 

</ 'J 



(27) 



Eqs. (I25ll-(l27b are inserted in Eq. (l24l l. giving Eyij which 
has three contributions (Evij)r', {Evij),^', and (Evij)0', corre- 
sponding to the three gradients of angular velocity in spherical 
coordinates. The total rate of energy dissipation in the surface 
layer is 



(28) 



'J 



In our computations the major contribution to E comes from 

(Evij)r'- 

3. Energy dissipation rates 

3.1. Sample E{(p',e') calculation 

We can compute Ejj, the rate of energy dissipation in the ele- 
ment j in a given parallel /, as a function of the azimuthal angle 
ip' of the element. Also, it is of interest to compute E,, the total 
rate of energy dissipation in the parallel /, as a function of the 
corresponding polar angle ff. As an example, in this section we 
show results from our model in the particular polytropic case 
with n=1.5 and using a surface layer depth ARi/Ri = 0.03, with 
corresponding average mass density p = 8.66x10"'' gr cm""*. 
The other input pai-ameters are Mi-5 Mq, M2-4 Mq, R\-3.2 
Rq, V = 0.003 Rl day an orbital period f = 10 days (a = 
40.62 Rq), and with two different values ySo = 1.2 and 2. 

Figure [1] shows the run of Ejj as a function of ip' in some 
selected parallels for Pq- 1.2. The upper curve corresponds to 
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log (sinS') 

Fig. 2. Rate of energy dissipation in the parallel / as a function 
of its polar angle, corresponding to the computations in Figure 
1 . Results with /3o - 1.2 and 2 are shown with empty and filled 
squares, respectively. The continuous lines show functions with 
a sin^ 0' dependence. 



the stellar equator, 0' - 90°, and following curves from top to 
bottom give results for the parallels with 0' ^ 81°, 73°, 64°, 
55°, 46°, 38°, 29° and 20°. A very similar results is obtained 
for the y0o=2 calculation, although the energy dissipation rates 
are significantly larger. 

Figure |2] shows the total rate of energy dissipation in the 
parallel /, as a function of its polar angle. Both cases /3o 
= 1.2 and 2 are shown in this figure, with empty and filled 
squares, respectively. The continuous lines show functions with 
a sin^ 0' dependence. Thus, our computations give the approx- 
imate dependence E, ~ sivJ 0', excluding polar angles in a 
region of approximately 10° around the poles. This implies 
that in our model (^)i has a sin^0| dependence, thus the 
dominant radial gradient in Eq. (l24t integrated in the volume 
of the parallel, which has a sin0^' dependence, gives £, ~ 
sin^ 0'- . Concerning this result, there is a theoretical study by 
Scharlemann (1981), which gives the tidal velocity field in a 
differentially rotating convective envelope of a component in a 
binary system in circular relative orbit, in the limit ySo — > 1 . For 
the case in which the star rotates uniformly, equation (40c) of 
Scharlemann (1981) gives the azimuthal component of the tidal 
velocity field. This equation implies, in our notation, (^), ~ 
sin^ 0', which diff'ers from our result. Yet we find in our compu- 
tations with /5q - 1.2, 2.0 that the surface layer does not rotate 
uniformly, contrary to the assumption in Scharlemann's equa- 
tion. Thus, the difference between the two results most likely 
stems from Scharlemann's assumption of uniform rotation over 
the stellar surface which, for /? 9^ 1 , is incorrect. 



3.2. Dependence of E on orbital separation 

In this section, we analyze the behavior of the total rate of en- 
ergy dissipation, E for binary system models with diff'erent or- 
bital separations and a variety of orbital eccentricities. 

3.2.1 . Circular orbits 

As in Toledano et al. (2007), we adopt as the test binary sys- 
tem one with masses mi = 5 Mq and TO2 = 4 M©; nii with a 
radius Ri - 3.2 Rq. The analysis was made for circular rela- 
tive orbits with several different values of orbital separation a. 
The energy dissipation rate given by Eq. (l2Ft was computed in 
OTi; we used 20 parallels distributed between the equator and 
polar angle 85°, with more than 10"^ surface elements, and a 
polytropic index n=1.5. In the coefficient of dynamic viscosity 
rj - vp the same approach as that of Toledano et al. (2007) was 
adopted. That is, the average mass density for the surface layer, 
ARi/Ri, is taken from a BEC stellar structure model compu- 
tation0 For this paper we used a modej| kindly provided by 
I. Brott (private communication, 2010). The value of the kine- 
matic viscosity v was taken as the lowest value allowed by the 
code (see discussion in Harrington et al. (2009) regarding this 
parameter). For the shortest period binary models with circu- 
lar orbits, this value is v = 0.003 /?| day \ which we adopted 
for the whole set of circular orbit models discussed in this sec- 
tion. Three different values of the thickness of the surface layer 
in mi were considered, ARi/Ri = 0.01, 0.03, 0.06|lThe cor- 
responding average mass densities are, 1.19xl0"^ 1.68x10"^ 
and 8.66x10"^ gr cm"^, respectively. The computations were 
made for two values /So - 1.2, 2.0 of the synchronicity param- 
eter. 

Figure|3]shows the energy dissipation E in nii as a function 
of the radius of the circular orbit. Empty symbols correspond 
to /?() = 1.2 and filled symbols to /5o - 2.0. Triangles, squares, 
and pentagons show results for ARi/Ri = 0.01, 0.03, and 0.06, 
respectively. A clear linear behavior is obtained on this log-log 
plot. For comparison, the continuous line in this figure shows 
a function with an a dependence on the orbital radius a. The 
data in the three cases ARi/Ri would shift vertically for differ- 
ent values of the mass density and viscosity, but the functional 
dependence on orbital separation remains the same. 

3.2.2. Elliptic orbits 

Elliptical orbit binary systems differ from those in circular or- 
bits in that the synchronicity parameter,/?, does not remain con- 
stant over the orbital cycle. This is because of the variation of 
the orbital angular velocity, Q. The computation takes into ac- 
count the changing value of /3, but we characterize each model 
by its synchronicity parameter at periastron, Pq. 

The energy dissipation rates for a grid of elliptical orbit bi- 
nary systems were computed following a similar procedure as 

The Binary Evolutionary Code, Langer (1991) 
mi= 5Mo, with an age of 4.29125e7 yrs, at which time its radius 
is 3.156 Rq. 

^ Appendix 1 shows the manner in which the E results depend on 
the chosen depth of the surface layer. 
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log (c/Ro) 

Fig. 3. Energy dissipation in mi as a function of the separa- 
tion a of the circular orbit. A polytropic index n = 1.5 is used. 
Empty symbols correspond to /3o -1.2 and filled symbols to 
jSo - 2.0. The correspondence with the layer thickness ARi/Ri 
is, in triangles: ARi/Ri - 0.01, in squares: ARi/R[ - 0.03, in 
pentagons: ARi/Ri = 0.06. The continuous line is a function 
with an dependence on a. 




e=0.0 a. 



^. 0.1 



1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 
log(a/Ro) 

Fig. 4. Energy dissipation rates (in ergs s"') from our model 
calculations for jSo -1-2 and different eccentricities as indi- 
cated. The dashed lines show the absolute value of E„rh given 
by Eq. A31 of Hut (1981), as written in our Eq. (gO). The dot- 
ted lines connect the points corresponding to each eccentricity. 
The short dot-dashed lines show the values of Ehm for e=0.8 
and lag angles q'=10° and 50°. 



in Section 13.2.11 for the mi + m2 - 5+4 M© binary system, 
but in this case eccentricities e =0.1, 0.3, 0.5 0.7 and 0.8 were 
assigned, and 20 latitude grid points were used in each hemi- 
sphere of mi. The remaining parameters were held constant: 
Ri =3.2, /3o =1.2, ARi/Ri =0.06, n =1.5, v =0.005 Rl/day. 
The changing separation of the two stars over the orbital cycle 
leads to a strong dependence of the energy dissipation rate on 
the orbital phase. Thus, in order to assign a value of E to each 
model calculation, an average value of the energy dissipation 
rate over orbital cycle, Eave, was computed by numerically in- 
tegrating E obtained at many phases distributed over the orbital 
cycle, and then dividing by the corresponding orbital period. 
The values of Eave are plotted in Figure |4] as a function of the 
major semi-axis of the orbit. 

Two aspects of this figure stand out: 1) there is a "family" 
of curves, one for each value of the eccentricity; 2) in each 
case, the scaling of Eave with a does not follow a unique linear 
(in the log-log plane) plot, as in the circular orbit case. Instead, 
for very short orbital periods, the slope is steeper and for very 
long periods, the slope is flatter than the Eg^e ~ relation for 
circular orbits. It is also interesting to note that the departure 
from this relation grows with increasing orbital eccentricity. 

3.3. Comparison with the "weak friction" modei 

Hut (1981) provided a derivation of the difi'erential equations 
for the tidal evolution under the "weak friction, equilibrium 
tide" model representation for the tidal interaction (Darwin 
1880; Alexander 1973). The fundamental assumption of this 
model is that the tidal bulges that are raised by the external per- 
turbing potential are only slightly misaligned with respect to 
the axis that connects the two stars. The system is then driven 



toward an equilibrium configuration through the action of the 
torques acting on the "retarded body". The misalignment of the 
bulges is caused by dissipation (i.e., friction) in the perturbed 
star. Hut's formalism includes high order terms in the orbital 
eccentricity, as required to properly assess the energy dissipa- 
tion rates in very eccentric systems]^ In this section we test our 
model by computing the energy dissipation rates for a 5-1-4 M© 
polytropic (n=1.5) binary system with a primary star's radius 
R\-3.2 Rq for a range in orbital separations and eccentricities. 
Each binary system of the grid is characterized by the energy 
dissipation rates averaged over the orbital cycle, Eave- We show 
that for a variety of orbital eccentricities, the behavior of 
is consistent with the predicted energy dissipation rates in the 
"weak friction" approximation given by Hut (1981) within the 
domain of applicability of this approximation. Our model also 
adequately predicts the pseudo-synchronization angular veloc- 
ity obtained by Hut (1981) for moderate eccentricities (e <0.3). 

3.3.1. Energy dissipation rates 

The rate of energy dissipation within the primary star associ- 
ated with the decrease in the orbital and rotation energy in the 
system is given by Hut (1981) in his Eq. A31 



■'Hut 



-3— G(mi +m2) 



mi, 



(R 



(1 -e^)-^^'^[A-2Bx + Cx^] 



^ The model has recently been revisited by Leconte et al. (2010) 
in the context of tidal heating in exoplanets. These authors extended 
Hut's equations to the case in which the orbital plane and the rotational 
equatorial plane of the objects do not coincide. A similar extension 
was also given in Eggleton et al. (1998). 



8 



Moreno, Koenigsberger & Harrington: Tidal flows and periastron events 



where A: is a constant that depends on the stellar structure, and 
X - w/n, with (JL> the angular rotation velocity of m\ and n - 
ripil - e^)^^^/(l + ef- the mean orbital angular velocity, with Up 
the orbital angular velocity at periastron. T - F?JGm\r, where 
T is the lag time, and 



1 + — e + 



255 



185 



25 



+ e + 

16 64 



, 15 2 45 4 5 6 



16 



C 



1 + + 4/. 



(29) 



After converting masses and distances to solar units, and writ- 
ing (jjjn in terms of rip, and recalling that ujup - j6o, the above 
equation can be rewritten as 



= -4.4 X lO'^^fer 



nil + m2\l m2 \ / "i i I o 



Mo l\Mel \Ro 
<(1 - eV^/^ [a - 2Bx + Cx^] , 



(30) 



with Ehui given in ergs s"' and x - /3o{l + e)^. 

For the lag time. Hut (1981; Eq. 4) gives t = a/ico - Q), 
where at and Q are the rotational and instantaneous orbital an- 
gular velocity, respectively, a is the lag angle, and is measured 
from the line that joins the centers of mi and 1112- This ex- 
pression may be re-written in terms of the instantaneous syn- 
chronicity parameter, /?, as t = a/Q.(J3 - 1). 

In eccentric binary systems, a is generally very small 
around periastron, but it can be quite large at other orbital 
phases|3 At the same time, the amplitude of the tidal bulges 
is largest around periastron and, given that this is when mi 
and m2 have their closest approach, the tidal interaction ef- 
fects are strongest. Hence, for the purpose of these exploratory 
calculations, we now make the following approximation: To 
compute the energy dissipation rates using Eq. (|30] |, we adopt 
T ~ aper/0.o(J3o - 1), whcrc ape,- is the lag angle at perias- 
tron which we assume is relatively small, and here tentatively 
adopt a=l°. We thus also use /3 = /3o- Using Kepler's relation 
Qo = G'/2(toi + m2yi^a-^^\\ + ef{\ - e2)-3/2 and converting 
to solar units. 



..i.(,xio-i4"f^r""-'>f (31) 



with Uper given in radians and r is in seconds. Note that when 
Eq. OTI ) is introduced into Eq. (|30] |. the dependence of Ehui on 
orbital separation ~ a thus flattening the slope from the 

relation that holds for a circular orbit. 

For the value of k, we use the expression given in Lecar et 
al. (1976; their eq. (2)) for the mass in the tidal bulge, M,: 



1 h 
M, = -kMi — 

2 Ri 



(32) 



An example of the orbital-phase dependence of a computed with 
our model may be found in Koenigsberger et al. 2003. 



where h is the height of the tide. Consider, for example, the 
model binary system with e=0.8, f =15 d, from the grid of mod- 
els discussed in the previous section. The height of the tide at 
periastron is ~0.03 R© (see below). The mass in the bulge is 
M, - f,Miayer, whcrc /, is the fraction of the mass in the outer 
layer, Miayer, that goes into to the bulge. For the particular case 
in question, Miayer = A-nR^ARi < p >=3. 67x10"'' Mq, where 
A/?i=0.06 Rq is the layer thickness and < p >=8. 66x10"^ g 
cm"^, is the average density of the layer. The maximum in the 
primary bulge approximately extends over 20° in azimuth, so 
// -10 Hence, this estimate yields k ~10"^. 

With these considerations, we use Eq. ( l30l l to estimate the 
energy dissipation rates in the "weak friction" approximation, 
Ehui, for the binary model parameters mi =5 Mq, m2-4 M©, 
Ri^3.2 R0,yS()=1.2, and k^SxlO'^. The dashed lines in Figure 
|4]correspond to the absolute values of Ehui given by Eq. (|30] |. 
The coincidence between the general behavior of these lines 
and that of the numerical computations is encouraging. Note 
that no relative shift or scaling is introduced, and yet the rela- 
tive separation between the curves corresponding to the differ- 
ent eccentricities is reproduced by our models. In addition, our 
model reproduces the general trend in the slope of the E^ve vs. 
a relation, particularly for large orbital separations. Thus, we 
find that our numerical calculation captures to a large extent 
the general scaling in the energy dissipation rates predicted an- 
alytically within the "weak friction" tidal model. 

However, although the scaling with e predicted by Eq. dSOl l 
is adequately captured by our model, the dependence on a ", 
with «=constant is not satisfied. Specifically, there is a marked 
excess in Eave for smaller orbital separations, with respect to 
Ehm- Several factors may be responsible for this difference. 
First, the assumption of a constant lag angle, a, is quite in- 
adequate. For example, assuming that for the e=0.8 calcula- 
tion a-lQ° and 50°, instead of the 1° used to draw the dash 
lines, significantly higher Eh„i values are obtained]^ This is il- 
lustrated by the short dashed lines in Figure |4] Second, the ap- 
plicability of the standard "equilibrium tide" model for high 
eccentricities is questionable (Efroimsky & Williams 2009). 
Finally, we note that our model assumes a spherical rigidly- 
rotating interior region, and that this approximation is less ap- 
plicable for small orbital separations, where the deformation of 
the star is greater. 

It is interesting to note that if E is plotted as a function of 
rper, the orbital separation at periastron, instead of the semi- 
major axis, the curves corresponding to the different eccen- 
tricities "collapse" toward a single curve. This is illustrated in 
Figure |5] Thus, the periastron distance is a more convenient 
parameter for the purpose of analyzing the average energy dis- 
sipation rates than the semi-major axis. 

3.3.2. Pseudo-synchronization angular velocity 

Pseudo-synchronization in an eccentric binary is a term used to 
describe the state that is closest to an equilibrium configuration. 
Hut (1981) demonstrated that the rotational angular velocity 



' We note that at periastron, a is indeed very small (see, for exam- 
ple, Figure [Tot, but grows after periastron passage. 
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Fig. 5. Energy dissipation rates of the previous figure but here 
plotted as a function of the periastron distance, r^^r- In this rep- 
resentation, the different curves for each eccentricity "collapse" 
onto a common region. The dashed lines show the absolute 
value of Eorb given by Eq. A31 of Hut (1981) for e=0.1 and 
0.8, also plotted as a function of rp^r- Symbols correspond to 
the different eccentricities: e=0.00 (open triangles), 0.01 (filled 
triangles), 0.10 (small pentagons), 0.30 (large pentagons), 0.50 
(crosses), 0.70 (stars) and 0.8 (filled pentagons). 




"Q°yooQQO 

05 1 15 " 0,5 1 1,5 

u/a, u/a. 

Fig. 6. Energy dissipation rates E„ye/10^^ < p > as a function 
of - (jj/Qo for four different eccentricities, e =0. 1 (top left), 
0.3 (top right), 0.5 (bottom left) and 0.8 (bottom right). The 
first three cases are for an orbital period P=6 d and the last is 
for 20 d. The ordinate is in units of 10^''' < p >ergs s"', where 
< p > is the average density. The dotted line indicates the level 
of minimum Eave which, for e=0.1 and 0.3 coincides with the 
analytical results of Hut (1981). 



for pseudo-synchronization, a>ps < Qo, where Qq is the orbital 
angular velocity at periastron. Specifically, 0.8 < Wp^/Qo < 10 
We computed sets of models with e=0.1, 0.3, 0.5 {P-6 d) 
and 0.8 (P-2Q d) holding all other parameters fixed except for 
/3o, which was varied from ~0.4 to the highest value admitted 
by the calculation[3Figure|6]illustrates the dependence of Eave 
on j6(). In all cases, there is a minimum in at a particular 
value (or range of values) of /3o. For e=0.1 and 0.3, the min- 
imum occurs around /3() ~0. 8-0.9, consistent with the "weak 
friction" approximation. For e=0.5 and 0.8, however, the min- 
imum occurs for 1.0< /3o <1.4, a result that is likely related to 
the caveats already mentioned at the end of the previous sec- 
tion. 



4. Synchronization timescales 

In this section we analyze the synchronization time scales, 
which will be defined below, using the model presented in 
Section |2] the procedure to compute the energy dissipation in 
Section 1231 and the values described in the previous sec- 
tion. The results will be compared with the earlier results of 
Toledano et al. (2007) and with a theoretical result given by 
Zahn (1977, 2008). 



Recall that in circular orbits, equilibrium configuration is at 
fio=a>/no=l, where £=0. 

' ' For the fixed value of v chosen for these calculations, very high 
stellar rotation rates cause neighboring surface elements to overlap, a 



4.1. Circular orbits 

As in Toledano et al. (2007), we now procede to estimate the 
synchronization time scale between the stellar rotation period 
and the orbital period. The synchronization time scale is the 
time needed by the system to arrive at a state in which both 
periods are equal, and this state is here assumed to be attained 
as a consequence of the energy dissipation in the surface layer. 
In our computations, at time f = the whole star, including 
the surface layer, is rotating with an angular velocity - 
/3o£Iq, with Qo the reference orbital angular velocity. For cir- 
cular orbits, Qo is directly the angular velocity in the circular 
orbit while for e >0 it is the angular velocity atperiastron. We 
adopt here as the condition for synchronicit}i[3 ^So = 1- Thus, 
with ySo > 1 there is an initial excess of rotational kinetic en- 
ergy over that in the synchronized state. As a first approxima- 
tion, ignoring the orbital change suffered by the binary system 
to reach the y6o = 1 state, the excess of rotational kinetic energy 
is proportional to (J3^ - l)Qo' "^^^^ the proportionality constant 
being 1/2, with / the moment of inertia. Thus for star mi the 
synchronization time is on the order of 



2E 



(33) 



From Kepler's law, £1^ has an a ^ dependence on the major 
semi-axis a. Thus, for a fixed value of y6o and with the approx- 



condition which halts the computation. 



Note, however, as shown by Hut (1981), the lowest energy state 
for the eccentric system corresponds to the pseudosynchronization pe- 
riod which, for moderate eccentricities, is /3 ~0.8-0.9. 
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Fig. 7. Synchronization time computed with Eq. (l33T l using the 
results given in Fig.|3] The meaning of symbols stays the same 
as in Fig. [3] The continuous lines are the theoretical result ob- 
tained by using Eq. (2.5) from Zahn (2008), under the assump- 
tion of a uniform density structure (/ = 2miR\l5), and with 
< V >=3.4xl0'^ cm^ s.'\ 



The top line is fory6o=1.2 and the 



bottom one for y6=2.0. 



imate a"^ dependence of E from Figure |3] Tsy„ will have an cfi 
dependence on a for circular orbits. 

On the other hand, Zahn (1977) has analyzed theoretically 
the tidal friction in close binary stars, also arriving at an 
dependence for stars in which a turbulent viscosity may be the 
main mechanism that produces the tidal friction. In a recent 
review of the topic (Zahn 2008), this expression for circular 
orbits is given as 



1 



'^Zahn ^dts 



1 (w-Q) ^niyR 

q 



(34) 



with q = m2/ini, I the moment of inertia of mi, O is the orbital 
angular velocity which for circular orbits is constant, t^iss = 
(oj - Q)/?j/aGmi is the dissipation time, with a the lag angle 
and G the universal constant of gravitation. Zahn (2008) sug- 
gests an approximation tjiss ~ R^l < v >, where < v > is an 
average kinematical viscosity, and we rewrite Eq. (l34l as: 



_ 5 X 10-^^ 



A) 



<v> \j6o - 1/ mi/Mo \/;, 



1 



Ki \ [a 



(35) 



where Tzo/m is in seconds. 

In order to compare our synchronization times obtained 
from Eq. (l33T l with those predicted in Eq. (l35T l, we now make 
the simplifying assumption of uniform density, from which 
/ = 2miR\l5. 

In Figure Q we show the values of t,,.,, obtained from our 
E values and Eq. (l33l l. These results have approximately the 
predicted dependence, as obtained in Toledano et al. (2007) 
with our earlier model. We also see that the derived synchro- 
nization times depend on stellar rotation, with the faster rotat- 
ing stars (J3q-2) attaining synchronization significantly faster 
than slower rotators (J3=l.2). 



The two continuous lines are Zahn's theoretical result given 
by Eq. i35[ for the two values of /3o, using nii - 5 Mq, 1112 - 4 
Mq, Ri - 3.2Rq, and < v >=3.4xl0'^ cm^ s"'. Our results are 
consistent with this theoretical result, except that the value we 
used for the kinematical viscosity in the calculations is signifi- 
cantly higher, v =0.003 d"' = 1.7xl0''* cm- s"'. It is impor- 
tant to note that in our calculation the tidal friction comes from 
the shear at the interface of each surface element with its neigh- 
bors and with the inner (rigidly rotating) boundary. This means 
that the action of the kinematical viscosity is limited to these 
shearing surfaces. In addition, because ours is a single-layer 
approximation, there is no buoyancy, and therefore our viscos- 
ity values are not associated with convection, as are those of 
Zahn (1977, 2008). Hence, at this stage, we are only able to 
draw an equivalence between the viscosity parameter used in 
our model and that used in Zahn's formalism through the com- 
parison of the synchronization timescales, as shown in Fig. [T] 
That < V >~ lO^^v is likely because in our calculations the 
whole energy dissipation is forced to occur in the surface layer, 
whereas in Zahn's formalism deeper layers are also involved in 
these processes. 

4.2. Eccentric orbits 

Eq. ( |33] | can be conveniently writterOl for an eccentric orbit in 
terms of the synchronicity parameter at periastron, /3o and the 
orbital eccentricity e, as 



TUyr) = 3.17 X 10 



38 



miR](/3l-m+e) 
Ea,eP(l-ey ' 



(36) 



where ;«i and Ri are the primary star's mass and radius, re- 
spectively, in solar units, and P is the orbital period in days. 
The synchronization times obtained with the above equation, 
using the Eat,,, values that are plotted in Figure |4] are shown in 
Figure |8] Several features of this plot are notable: 1) there is a 
clear difference between moderately eccentric (e <0.3) and the 
very eccentric systems; 2) at large orbital separations (a >100 
Rq), the slope significanlty flattens from the r,,,, ~ (a/Ro)^ 
that applies for circular orbits0; and 3) for a fixed value of a, 
the timescale for synchronization of the eccentric binaries is 
shorter than the corresponding circular orbit binary, with dif- 
ferences as large as five orders of magnitude visible in Figure 

m 

The separation into these different "families" is eliminated 
in part by using the periastron distance, rpe,-, instead of the 
semi-major axis of the orbit. Figure |9]illustrates the same val- 
ues of E shown in Figure [8] but in this case rj,.„ is plotted as a 
function of r^^r- The different curves "collapse" onto a common 
region in the logijsyn) vs. log{rper) plane, with significant dif- 
ferences occuring only for very large orbital separations. This 
is consistent with the fact that, although the mean distance be- 
tween the two stars increases with e, the distance at periastron 



" Using Oo = 1-428 x I R^Po, and v„„=y8oi?i(l + 

eyi^l[Q.Q2P{\ - efl\ where Vto, is in km s , P is in days, and 
inRo. 

Note that for a >100 Rq, t.„,„ ~ {ajRef -^^^, similar to the syn- 
chronization timescale predicted by Tassoul (1987). 
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Fig. 8. Synchronization from the energy dissipation rates com- 
puted with our model using a kinematical viscosity v=0.005 R| 
d ' and Eq. (|36] |. The different "families" are labeled by their 
corresponding eccentricity. The solid line is the synchroniza- 
tion timescale from Zahn's (2008) relation, assuming a con- 
stant density structure and < v >=1.4xl0'^ cm^ s"'. This 
relation is valid only for e ~0.0. The dotted lines illustrate 
the trends for higher eccentricities from our computed mod- 
els. Note that Tsv„ ~ a*-^^^ for high value of e and a. is in 
units of years. 
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Fig. 9. Same synchronization times in the previous figure, but 
here plotted as a function of periastron distance, Vpgr- The dot- 
ted lines connect the points corresponding to each of the eccen- 
tricities. TjyB is in units of years. 



strongly decreases, and as emphasized by Leconte et al. (2010), 
most of the work caused by tidal forces occurs at this point of 
the orbit. 




Fig. 10. Top: maps of radius residuals (r((^')-/?i) for e=0.8, 
P=15d,yS() =1.2 at periastron (left) and 0.833 days after (right). 
The sub-binary longitude ((^' =0) is at the far left on each map 
and (/)' = 180° is at the center. The sense of the stellar rotation is 
from left-to-right. White/black represents maximum/minimum 
residuals. Bottom: corresponding plots of the radius residuals 
along the equator. 

The actual values of t,,.,, depend on the values of v and k. 
But the scaling of t,,.,, with separation at periastron is as plotted 
in Figure[8] unless v and k depend on orbital separation!'^! 

5. Orbital phase-dependence and periastron 
events 

Consider the binary system with e=0.8 and P=15d from the 
grid of models described in the previous sections. With /3()-1.2, 
the primary star, mi, is rotating at 193 km s"'. The maps in the 
top panels of Figure 10 are a color-coded Mollweide projection 
of the stellar surface, with white corresponding to the maxi- 
mum elevation from the unperturbed stellar surface and black 
to the largest depression below the unperturbed stellar surface. 
Azimuth angle (p-Q°, corresponding to the location of m2, is on 
the left edge of the map, i^=360° on the right edge, and ip-lSO° 
is at the center. At periastron (Figure 10, left), the surface ac- 
quires a shape that is typical of that assumed in the "equilib- 
rium tide" approximation. That is, a tidal bulge is raised on both 
hemispheres with the larger bulge pointing in the general direc- 
tion of the companion, m2, but "lagging" slightly behind. The 
maximum height of the primary bulge at the equator (Figure 
10, bottom left) is -0.03 RqEI 

The "equilibrium tide" configuration does not subsist very 
long, however. The right-hand panel of Figure [TO] illustrates 
the deformations just 0.833 day after periastron, and clearly 
the "equilibrium tide" approximation no longer represents the 

Because v is generally associated with turbulence, and because 
the perturbations of the surface velocity field in short-period binaries 
is very strong, the turbulent component in these systems may be much 
stronger than in wider binary systems. Thus, the use of a constant 
viscosity for all orbital separations may not be the most appropriate. 

A rough estimate (c.f., Zahn 2008) of the tidal bulge eleva- 
tion for a non-rotating star with similar characteristics, SR/Ri ~ 
{milm2){Rilrp„)* ~ 0.02, thus predicting SR=Q.Q1 Rq, a factor 2.5 
larger than predicted by our model. 
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Fig. 11. Top: maps of tangential velocity residuals for e=0.8, 
P=15d,j8o =1.2 at periastron (left) and 0.833 days after (right). 
The sub-binary longitude (^' =0) is at the far left on each map 
and (^' = 180° is at the center. The sense of the stellar rotation 
is from left-to-right. White/black represents positive/negative 
residuals; i.e., particles moving faster/slower than the unper- 
turbed stellar rotation velocity. Bottom: corresponding plots 
along the equator. 

stellar surface shape. Now it is best described in terms of a 
large number of smaller-scale spatial structures. These would 
fall into the "dynamical tide" representation, if they were as- 
sociated with the normal modes of oscillation of the stellar 
interiorFl The amplitude of these structures declines with or- 
bital phase as apastron is approached. A more in-depth analy- 
sis of the smaller spatial structures and their dependence on the 
model input parameters goes beyond the scope of the current 
paper, but the important point to note is that the "equilibrium 
tide" approximation is valid only around periastron passage for 
this particular highly eccentric model binary case. We find that 
for systems in our computed grid with low eccentricities and/or 
large orbital periods, the approximation remains valid over the 
orbital cycle, but not in the more extreme cases. However, even 
there, the tidal bulges raised by the secondary cannot be as- 
sumed to be steady over the orbital cycle, neither in amplitude 
nor in location. 

We illustrate the azimuthal velocity field, v^-, in Figure [TT] 
The explanation given by Tassoul (1987) in the context of cir- 
cular orbits is most appropriate to describe the phenomenon 
displayed at periastron (left): "Assume that at some initial in- 
stant the primary is rotating rigidly with a constant supers yn- 
chronous angular velocity. If there were no companion, the sys- 
tem would remain forever axisymmetric. Because of the pres- 
ence of the secondary, however, a fluid particle moving on the 
free surface of the primary will experience small accelerations 
and decelerations in the azimuthal direction. Hence, in the sec- 
tors Q° < (f' < 71 12 and n < if' < 37t/2, the azimuthal compo- 
nent v^, of this fluid particle will be smaller than the typical 
value by a small amount. On the contrary, in the other sec- 
tors, v^ will be larger" Note that in the case shown in Figure 
[TT]the pattern of positive and negative residuals is slightly ro- 
tated with respect to Tassoul's description because this is not a 

" Recall that the stellar interior here is a rotating, rigid body. 




Fig. 12. Top: maps of energy dissipation rate for e=0.8, P= 15d, 
ySo =1.2 at periastron (left) and 0.833 days after (right). The 
sub-binary longitude (^' =0) is at the far left on each map and 
(/?' = 180° is at the center The sense of the stellar rotation is 
from left-to-right. White/black represents maximum/minimum 
residuals. Bottom: corresponding plots along the equator. 

circular-orbit case, and also because of the "lag" of the primary 
tidal bulge. This quadrupolar pattern rapidly breaks down into 
smaller scale structures after periastron, corresponding to the 
smaller spatial structures of Figure [TO] (right). 

Highest energy dissipation rates are associated with the re- 
gions having large perturbations in the azimuthal motion, both 
positive and negative. At periastron (Fig. [12] left), the high- 
est E values occur close to the sub-binary longitude. But by 
0.833 day later (Fig. [12] right), the large regions where most 
of the E is generated have broken down into small spatial scale 
structures. It is noteworthy, however, that some of these smaller 
structures are associated with significantly high values of E. 
Also noteworthy is that nearly the entire E occurs at +45° from 
the equatorial latitude. 

The main conclusions we draw from this example are that 

1) periastron passage causes a sudden and brief change in the 
surface perturbations, leading to significant dynamical effects; 

2) the tidal shear energy dissipation rates are concentrated (in 
this example) around the equatorial latitude. Hence, the ques- 
tion arises as to whether these phenomena may provide a nat- 
ural explanation for the periastron events described in the in- 
troduction, as well as the creation of a disk from the ejected 
material. The topic of how the material may be ejected goes 
beyond the scope of this papeQ but a plausible scenario is 
that additional energy is fed into sub-photospheric layers by the 
tidal shear, combined with the possible additional turbulence 
caused by the rapidly-changing motions on the stellar surface. 
We suggest that one of the crucial ingredients is the extremely 
rapid and strong forcing to which the stellar surface is subj ected 
during periastron passage in very eccentric systems. Figure [T3] 
shows a plot of £ as a function of orbital phase for the e=0.8, 
P=15 d binmy model discussed above and two other cases for 
comparison]^ The rapid rise in the E value at periastron is ev- 
ident. 

'** See, for example, Lee (1993) 

" Note that actual binary systems with this eccentricity and orbital 
period are extremely rare or non-existent, as can be seen in Fig. 1 of 
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Fig. 13. Tidal shear energy dissipation rate for very eccentric 
binary systems as a function of orbital phase measured from 
periastron passage. E is in ergs s ' and p in gm cm"^. Shown in 
this figure are two cases with e=0.8 (P=15 and 350 days) and a 
case for 6 Sco. Note that peristron passage is always associated 
with a very large time variation of E. 

The second case depicted in Figure [T3] corresponds to a 
e=0.8, P=350 d model from our grid, and the unperturbed equa- 
torial rotation velocity of mi is ~8 km s"'. The E values un- 
dergo a ~6 order-of-magnitude variation between apastron and 
periastron, most of which occurs within +70 days of perias- 
tron passage. Finally, the third case illustrated in Figure [13] 
is an exploratory model calculation for the extremely eccen- 
tric (e—0.94) Be-type binary system 6 Scol^ The primary star 
in this system has an estimated equatorial rotation velocity of 
~200 km s"', which means that it is rotating extremely fast, 
compared to the orbital angular velocity of m2. The star un- 
dergoes an increase in E by ~9 orders of magnitude between 
apastron and periastron, with a very steep rise to maximumjust 
15 days before periastron. 

6. Discussion 

A large body of research, starting with Darwin (1880), exists on 
the effects produced by tidal interactions on the long-term evo- 
lution of binary systems. The primary concern of these inves- 
tigations is generally the timescales for circularization of the 
orbit and synchronization of the rotational and orbital angular 
velocities. The three general approaches that have been devel- 
oped for analyzing the effect of tidal dissipation in these inves- 
tigations are 1) the "equilibrium tide" model; 2) the "dynami- 
cal tide" model; and 3) the hydrodynamical model of Tassoul 
(1987). 



Mazeh (2008), which is probably because of the very large E, which 
should rapidly tend to circularize the orbit. 

The parameters used for this calculation are mi= 15 M©, m2= 2 
Mo, P=3944.7 d,Ri=7 Ro,ySo=28.5, ARi/Ri=0.06, /3o=2S.5, v=0.08 



In the "equilibrium tide" model (Darwin 1880; Alexander 
1973; Hut 1980, 1981; and see Eggleton, Kiseleva & Hut 1998 
and references thereinj^ the shape of the stellar surface is ap- 
proximated to that which it would adopt in equilibrium under 
the effect of the gravitational-centrifugal potential of the sys- 
tem. The viscous nature of the stellar material leads to a mis- 
alignment (refereed to as "lag") between the axis that connects 
the two tidal bulges and the line that connects the centers of the 
two stars. The action of the external gravitational field on this 
misaligned body leads to a torque which, in turn, leads to an 
interchange of energy and angular momentum between stellar 
rotation and orbital motion. Dissipation by viscosity results in 
the decrease in the eccentricity of the orbit and the synchro- 
nization of the rotational and orbital motions. It is important 
to note that the "equilibrium tide" approach implicitly assumes 
that the star goes through a succession of rigidly rotating states, 
until the equilibrium configuration is reached. 

The "dynamical tide" approach is based on the assump- 
tion that the star may be viewed as an oscillator whose nor- 
mal modes may be excited by the variation in the gravitational 
potential. Among the first mathematical formulations available 
in the published literature are those of Fabian, Pringle & Rees 
(1975), Zahn (1975), and Press & Teukolsky (1977). In more 
recent years, a large amount of effort has been invested in fur- 
ther developing these concepts, as reviewed by Savonije (2008) 
(see, also McMillan, McDemott & Taam 1987; Goldreich & 
Nicholson 1989; Dolginov & Smel'Chakova 1992; Lai 1996; 
Mardling 1995; Kumar & Goodman 1996; Kumar & Quataert 
1997). Three types of modes are believed to be excited a) 
the internal gravity modes, where the restoring force is the 
bouyancy force in stably stratified regions (Zahn 1975; Zahn, 
1977; Goldreich & Nicholson 1989; Rocca 1989); b) the acous- 
tic modes, where the restoring force is the compressibility of 
the gas; and c) the inertial modes where the restoring force 
is the Coriohs force (Witte & Savonije, 1999a, 1999b). These 
modes are the mechanism by which energy from the orbital 
motion is absorbed by the star. The system is driven toward 
the equilibrium state by the damping of these modes. Hence, 
it is necessary to make an assumption concerning the damping 
mechanism. In stars with an outer convective layer, it is gen- 
erally assumed that turbulent viscosity is the dominant mecha- 
nism. In stars with outer radiative layers, radiation damping is 
assumed to be the dissipative mechanism. It is interesting that 
in an eccentric system, the gravitational forcing goes through 
a large range of frequencies as the orbit progresses from apas- 
tron to periastron. Hence, there is a good chance that one of 
these frequencies will resonate with the star's normal modes, 
thus exciting the particular oscillation mode in the star. At the 
same time, as periastron is approached, however, the star be- 
comes more deformed, thus potentially modifying its normal 
modes. The above methods rely on a mathematical formulation 
of the star's response to the perturbing potential of its compan- 
ion, and are usually treated in the linear regime. Thus, they are 
generally limited to slow rotational velocities and/or low orbital 
eccentricites. 



^' An illustrative description and discussion of the limitations of this 
model may be found in Ferraz-Mello et al. (2008). 



14 



Moreno, Koenigsberger & Harrington: Tidal flows and periastron events 



The third approach, introduced by Tassoul (1987) for early- 
type binaries, involves large-scale hydrodynamical motions 
within the tidally distorted star. It is based on the assumed pres- 
ence of large-scale meriodional circulation, superposed on the 
mean azimuthal motion of stellar rotation. This circulation is 
believed to arise because a fluid particle moving on the free sur- 
face of the primary will experience small accelerations and de- 
celerations in the azimuthal direction, depending on its location 
with regard to the sub-binary longitude|3 These meridional cir- 
culation currents transport angular momentum very efficiently, 
which leads to very rapid synchronization timescales; i.e., the 
spin-down time primarily depends on (a/R)^ '^''. This mecha- 
nism is thus a more long-range mechanism compared with the 
other two. The method recognizes that near the surface, the 
Coriolis force and the viscous force must play a donominant 
role in the mechanical balance. Hence, this approach signifi- 
cantly differs from the other two in that it captures the hydro- 
dynamic processes near the stellar surface. 

Contrary to the three methods described above, where the 
primary concern is the long-term (>10^ yr) evolution of the 
orbital and stellar rotation parameters, our method was devel- 
oped to analyze the tidal effects on orbital timescales (~days). 
However, we are able to compare the results of our model with 
results from the other models through the energy dissipation 
rates which, in turn, allow an estimate to be made of the syn- 
chronization timescales. We find that the results of our model 
for circular orbits are fully consistent with the results of Zahn 
(1977, 2008) for stars with convective envelopes. We also find 
that for eccentric orbits, our results for the scaling of the en- 
ergy dissipation rates with orbital eccentricity and separation 
are consistent with the predictions of Hut (1981). 

Zahn's synchronization timescale for convective envelopes 
is based on the "dynamical tide" formalism, while Hut's 
method is based on the "equilibrium tide" approximation. Our 
method does not rely on any a priori assumption regarding the 
mathematical formulation of the tidal flow structure. Instead, 
we solve the equations of motion of a thin, deformable sur- 
face layer that is allowed to respond to the gravitational field 
of mi and m2, the Coriolis force, the centrifugal force, the gas 
pressure and the viscous forces from the shearing flows. The 
numerical solution of this system of equations yields the veloc- 
ity field on the stellar surface, from which the energy dissipa- 
tion rates may be computed. Thus, ours is a computation from 
first principles, the only simplifying assumptions being those 
which are inherent to the one-layer approximation. In particu- 
lar, our method neglects the deformations produced by the tidal 
efifects on the layers below the surface layer However, because 
the tidal forces decline as the tidal effects are weaker in 
deeper layers of the star than on the surface (see, for exam- 
ple, Dolginov & Smelch'akova 1992). On the other hand, we 
are fully aware that the detailed behavior of underlying layers 
will be coupled to that of the surface layer Thus, at this stage 
our results are merely indicative of the detailed behavior of the 
surface layer as a function of the orbital phase. 



Indeed, our simulations show that these azimuthal motions are 
the dominant cause of photospheric line-profile variability on orbital 
timescales (Harrington et al. 2009). 



Among the advantages of our method is that it allows a 
qualitative assessment of the response of the stellar surface 
under arbitrary conditions and it captures the essence of the 
three general approaches described above. That is, the solu- 
tion of the equations of motion yields a first-order deforma- 
tion attributable to the equilibrium tide and a lag angle; a sec- 
ond order deformation attributable to dynamical effects!^ and 
because of the inclusion of gas pressure, it captures some of 
the hydrodynamical effects. Interestingly, our results for the 
synchronization timescales in eccentric orbits are consistent 
with the equilibrium tide approach at intermediate orbital sep- 
arations, but suggest much faster synchronization timescales 
for larger orbital separations, consistent with Tassoul's hydro- 
dynamic model0 From the observational standpoint, Abt & 
Boonyarak (2004) found that B and A stars in binaries with 
periods as long as 500 days have rotational periods signifi- 
cantly shorter than the corresponding single stars. They there- 
fore concluded that the synchronization process can have an 
inpact on the stellar rotation even for relatively wide binaries 
with periods as long as 500 days, thus supporting the presence 
of a longer-range mechanism. Claret et al. (1995) and Claret 
& Cunha (1997) also have noted that the Zahn models and the 
Tassoul models each can explain the observational data only in 
part. 

The major limitation of a one-layer model lies in neglect- 
ing the oscillatory nature of the stellar interior. In particular, 
the interaction of the normal modes of oscillation of the star 
with the external forcing frequency provides a means by which 
orbital energy might be captured by the star and deposited as 
thermal energy. Hence, we cannot at this time probe interest- 
ing situations associated with resonances. On the other hand, 
in long-period, very eccentric binary systems the peristron pas- 
sage involves very large temporal gradients in the surface prop- 
erties of the star The tidal shear energy dissipation rates, in 
particular, change dramatically over a very small portion of the 
orbital cycle. Our exploratory calculation for the 6 Sco binary 
system (orbital period P=3944.7 days) indicates that E changes 
by five orders of magnitude over the 82 days before perias- 
tron, two of these occurring over the 15 days before peristron. 
The combination of added energy, extended radius, and high 
surface flow velocities may be envisaged to produce a destabi- 
lizing effect that leads to a mass-ejection event. Interestingly, 
we expect that the ejecta should form a disk-like structure, be- 
cause the strongest tidal effects occur in the equatorial plane. 
Observational evidence for such a disk-forming event actu- 
ally exists for the previous orbital cycle (Miroshnichenko et 
al. 2003). Confirming this behavior in the upcoming periastron 
passage of 201 1 would strengthen the hypothesis. 

Finally, we note that it may be productive to use the peri- 
astron distance, rpc,-, instead of the semi-major axis, a, when 



Clearly, the dynamical effects predicted by the model refer to the 
surface layer only. 

Note, however, that our model does not currently include merid- 
ional motions. The connection with the Tassoul model resides in the 
concept of the azimuthal flow velocities. These have significantly 
faster speeds than those of meridional flows and act on very short 
timescales. 
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comparing statistical studies of observational properties with 
theoretical models of tidal evolution. We show that by using 
rper, the synchronization timescales for binaries with different 
orbital eccentricities "collapse" into a small region in the t,,,,, 
vs Tper diagram, whereas on the r,,,! vs a diagram, there is a 
different "family" of curves corresponding to each eccentric- 
ity. This means that for a fixed value of a, the synchronization 
timescale derived from eccentric binary systems is significantly 
shorter (by orders of magnitude) than a corresponding e-Q sys- 
tem. However, if rper is used instead, the spread between sys- 
tems with different eccentricites is smaller. 

We emphasize moreover that there is a significant depen- 
dence of Tsyn on the stellar rotational velocity. However, it must 
be kept in mind that the important parameter is not only the 
equatorial rotational velocity, but also the synchronicity param- 
eter. Hence, slowly rotating stars may undergo strong tidal 
interaction effects if, for example, /?o << 1- Whether strong 
tidal effects can contribute to alter chemical abundances is a 
question that merits examination, because this could potentially 
be an explanation for some of the discrepancies recently found 
between stellar structure evolutionary models and observations 
(Hunter et al. 2008; Brott et al. 2010). 

Clearly, many of the outstanding questions associated with 
tidal interactions require a more complete treatment than that of 
any of the currently developed models. In particular, the depar- 
tures from spherical symmetry imply that a one-dimensional 
treatment of the radiative transfer processes is inadequate. 
Furthermore, the energy that is dissipated due to tidal interac- 
tions, regardless of the mechanism invoked, is an added source 
of energy that is deposited in sub-photospheric layers of the 
star and needs to be considered when analyzing energy trans- 
port processes and determining the structure of the outer stellar 
layers. 

7. Conclusions 

We derived a mathematical expression for computing the rate 
of dissipation, E, of kinetic energy by the viscous flows that 
are driven by tidal interactions and used it to explore the effects 
caused by periastron passage in eccentric binary systems. Our 
method does not rely on any a priori assumptions regarding 
the mathematical formulation for the perturbations in the star 
caused by its companion. We show that our results for the scal- 
ing of E with orbital separation for circular orbits are the same 
as those derived analytically by Hut (1981), and that our re- 
sults yield the same scaling of synchronization timescales with 
orbital separation as found by Zahn (1977, 2008) for convec- 
tive envelopes. For elliptical orbits, our results reproduce the 
dependence of E on the eccentricity and the general trend with 
orbital separation as predicted by Hut (1981). In addition, the 
model also adequately predicts the pseudo-synchronization an- 
gular velocity for moderate eccentricities (e <0.3), as obtained 
by Hut (1981). 

We explored the distribution of E over the stellar surface, 
and found that its configuration in small eccentricity binaries 
is consistent with the quadrupolar morphology of the "equilib- 
rium tide". That is, maximum E values tend to be located at 
four longitude regions where the fastest azimuthal motions are 



present. In very eccentric systems, however, the morphology 
significantly changes as a function of orbital phase. Although 
generally quadrupolar-like at periastron, shortly after this time 
E is distributed into smaller spatial structures, which suggests 
an analogy with the "dynamical tide" models. In addition, the 
perturbations are in general strongest around the equatorial lat- 
itude |3 

We speculate that the very large and sudden increase in E 
that occurs around periastron, combined with the rapid growth 
rate of the velocity perturbations is a promising mechanism 
for explaining the increased stellar activity and mass-ejection 
events often observed around this orbital phase. Two factors 
may be involved in causing the periastron events. First, E 
is an additional source of energy that is deposited in sub- 
photospheric layers and represents a time-dependent pertur- 
bation to the temperature and differential velocity structure of 
these layers. 

Second, in very eccentric binaries, the growth rate of the 
surface perturbations increases very abruptly at periastron. For 
example, in the very eccentric binary system 6 Sco (e=0.94, 
f =3944.7 d), our exploratory calculation indicates that E 
changes by five orders of magnitude over the 82 days before 
periastron, with a two-orders of magnitude increase over the 
15 days before peristron. The qualitative phenomenon is one 
in which portions of the stellar surface are forced to suddenly 
rotate faster or slower than their neighboring inner layers , thus 
leading to a variable differentially-rotating structure which sig- 
nificantly differs from its unperturbed structure. Furthermore, 
most of the activity occurs around the equatorial latitude. The 
consequence for the stability of the outer layers under these 
conditions is a problem that merits further attention. 
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Table A.l. Summary of the 5+4 Mq model binary system mod- 
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2-265 


1.2 




3 


0.1 


3-255 


1.2 




4 


0.1 


6 


0.4-2.1 


minimum £„,,<> at ySo=0.90 


5 


0.1 


6 


1.2 


V =0.001-0.088 R| day-' 


6 


0.1 


6 


1.2 


Ai^i/i?! =0.02-0. 10 


7 


0.3 


3-255 


1.2 




8 


0.3 


3-390 


2.0 




9 


0.3 


6 


1.2 


Ai?i//?i=0.02-0.10 


10 


0.3 


6 


0.4-2.35 


minimum E^^g at /Jo =0.80 


11 


0.5 


4-48 


1.2 




12 


0.5 


6 


0.4-2.35 


minimum at /3o= 1 .05 


13 


0.7 


8-255 


1.2 




14 


0.8 


13-350 


1.2 




15 


0.8 


20 


1.2 


ARi/Ri =0.02-0.10 


16 


0.8 


20 


0.4-1.80 


minimum Eave at/3o ~1.15 


17 


0.8 


30 


0.4-1.25 


minimum at A ~1.15 



Appendix A: Grid of model calculations 

Table 1 lists the input parameters of the grid of model calcula- 
tions used in this paper. The first column identifies the set num- 
ber; the second column lists the eccentricity; the third column 
lists the range in orbital periods included in the given set; the 
fourth column lists the value of /3o used in the calculation, or 
the range in /3q values; and the last column contains comments. 
Unless otherwise noted, all calculations were performed with 
a kinematical viscosity parameter, v=0.005 RQ/day=2.8xlO 
cm-/s. 
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Appendix B: Dependence of Eave on layer depth 

We examined the manner in which the energy dissipation rates 
depend on the choice of the depth {ARi/Ri) of the surface layer 
that is modeled. This was done holding all parameters fixed, 
except for the values of dRi/Ri, which were varied from 0.02 to 
0.10. The results for e=0.1, 0.3 and 0.8 are illustrated in Figure 
IB. II which show that for the models we have run for this paper, 
the calculations using different layer depths yield values of Eat.^ 
that differ by no more than a factor of two. 
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Fig.B.l. Dependence of on the choice of layer depth for 
models computed with e=0.1, 0.3 and e=0.8, and correspond- 
ing orbital periods as listed in the label. The average density 
for each layer depth used for the E^ve calculation were ob- 
tained from the stellar structure model of a 5 Mo, 3.15 Ro star 
(I. Brott, private communication, 2010). This shows that for a 
fixed set of binary parameters, the energy dissipation rates de- 
pend on the choice of layer depth only within a factor of ~2. 



